Wellbore Positioning System and Method

ABSTRACT

A method and a system are provided for determining the relative positions of a wellbore and an object. The wellbore is represented by a first ellipse and the object is represented by a second ellipse. The first ellipse represents the positional uncertainty of the wellbore and the second ellipse represents the positional uncertainty of the object. The method includes receiving input data relating to a measured or estimated position of the wellbore and the object. In addition, the method includes calculating an expansion factor representing an amount by which one, or both, of the first ellipse and the second ellipse can be expanded with respect to one or both of respective first and second sets of elliptical parameters so that the first and second ellipses osculate. Further, the method includes determining, based on the calculated expansion factor, position data indicative of the relative positions of the wellbore and the object.

TECHNICAL FIELD

The present invention relates to a computer-implemented method and asystem for determining the relative positions of a wellbore and anobject. Spatial relationships between two ellipses, each of whichrepresents the positional uncertainty of a wellbore, are utilized todetermine the conditions governing osculation between the two ellipses,expressing the determination as an expansion scale factor.

BACKGROUND

As the drilling of a wellbore in hydrocarbon reservoir (for example, anoil or gas reservoir) proceeds, the positional uncertainty at any pointin a well is dependent on a number of factors, including the positionaluncertainty of the surface location, the well's geographical locationand trajectory and the various instruments used to survey the well. Bypositional uncertainty is meant positional uncertainty of the well'sgeographical location, positional uncertainty of its trajectory etc. Theexpected behaviours of these instrument types are presented asinstrument performance models. Application of these models quantifiesthe uncertainty of the true wellbore position for a stated confidence.

Referring to FIG. 1, the positional uncertainty 1 about a pointrepresenting the calculated position of the centre of a wellbore iscommonly represented as an ellipsoid with its principal axes alignedwith the high-side, right-side and along-hole directions. In thiscontext, high-side is the direction normal to the wellbore in thevertical plane and right-side is the direction normal to both thewellbore direction and the high-side and so lies in the horizontalplane. The ellipsoid usually also accounts for the dimensions of thecasing or open hole of the wellbore. The size of the ellipsoid variesaccording to the wellbore trajectory shape, survey instruments used incalculating the position, and selected confidence limits.

Using this model, at any time and point in space, which is to saypositions down a wellbore and its surrounding volume, the resultingpositional uncertainty about a wellbore along its trajectory is theenvelope of the ellipsoids; a curved, continuous cone 2 with a curvedend. The interference between two adjacent wells can be visualised asthe interference between the two cones. The positional uncertainty canchange over time as data is re-processed or more data is acquired. Italso changes when the wellbore is resurveyed using a more accurateinstrument system, for example when a high accuracy gyroscope is run ata casing point 3. This then narrows the cone, as shown in FIG. 1.Therefore, if a new measurement is taken at a subsequent point along thetrajectory of the wellbore, the positional uncertainty decreases, thenincreases as the distance from the measurement point increases.

To a good approximation, at any given point along the wellbore theintersection of a plane that is normal to the along-hole direction ofthe wellbore with the cone can be represented as an ellipse. Therefore,the problem of calculating the interference between two wells can bereduced to that of calculating the distance between two such ellipses.This simple geometrical model has been adopted by various standardsorganisations to define minimum acceptable separation distances betweentwo wellbores, for example the Norwegian “Norsk SokkelsKonkurranseposisjon” (NORSOK) D-10.

The separation between wellbores in 3D can be represented in 2D using acollision avoidance plot, also known as a travelling cylinder or normalplane diagram. In this representation the intersection of, for example,an existing and a planned well (or two planned wells) is displayed on aplane, constructed normal to the planned well. The planned well is keptat the centre of the plot and therefore the relative separation betweenthe planned well and the adjacent well is indicated by the locus ofpoints obtained at successive depths. At any point in the subject wellthe plane also intersects the curved cone and at low or modest angles ofincidence between wells the intersections with the cones appears, togood approximation as two ellipses. During drilling the as-drilled andprojected positions are shown on the same plot. The planned well is alsoreferred to as the subject or reference well.

If x and y are orthogonal coordinates in the normal plane then theseparation δ between the wells can be calculated using Eqs. 1 to 3,below. In the absence of bias this is also the separation between theerror ellipses. Further adjustments can be made if required.

Δx ₀ =x _(0,2) −x _(0,1)  (1)

Δy ₀ =y _(0,2) −y _(0,1)  (2)

δ=(Δx ₀ ² +Δy ₀ ²)^(1/2)  (3)

In practical terms the minimum approach distance δ_(min) between thewellbores must be greater than the sum of the open hole and casingradii, δ_(min)>(d_(h)+d_(c))/2, where d_(h) is the hole diameter andd_(c) is the casing diameter. This criterion automatically satisfies themathematical constraint δ≠0.

Currently, the relationship between two adjacent ellipses isapproximated as a “separation factor”, k_(s). In this representation theellipses are related only by the line passing through their centres.Because of this, the calculation of the characteristic length s for eachellipse may be performed independently of the other. Two common methodsare the centre vector method (CVM) and pedal curve method (PCM).

$\begin{matrix}{k_{s} = \frac{\delta}{s_{1} + s_{2}}} & (4)\end{matrix}$

Because of mathematical difficulties, existing methods for calculatingseparation factors are approximations and may be either too optimisticor too conservative, particularly for ellipses with high eccentricities.

For example, FIG. 2 shows how the currently used “centre vector method”(CVM) is used to calculate a separation factor between two wellbores. Inthe CVM the characteristic lengths s₁ and s₂ are determined from thepoint of intersection of each ellipse, marked A and B in FIG. 2, withthe line δ joining their centres. The separation factor k_(CVM) iscalculated using Eq. 4. In this case the ellipses extend beyond theirpoints of intersection and will touch before the separation factorreaches unity. Therefore, separation factors calculated using thismethod may be too optimistic. Such overly optimistic calculations of theseparation factor can lead to safety issues when planning and drillingwells based on computer simulations.

FIG. 3 shows an alternative method of calculating a separation factorbetween two wellbores, the “pedal curve method” (PCM). In the PCM, thecharacteristic lengths s₁ and s₂ are determined from the line that isboth tangent to the ellipse and is orthogonal to the line δ joiningtheir centres. The first step is to determine the points of tangency,marked A and B in FIG. 3. In this case the tangent lines meet and theseparation factor k_(PCM) reaches unity before the ellipses touch.Therefore the separation factors calculated using this method may be tooconservative, leading to unnecessary shut-in of wells or missedopportunities.

Although the separation factors calculated by either the centre vectoror pedal curve methods are relatively easy to calculate, neither methodis a faithful representation of the geometrical relationship between thetwo ellipses. As shown in FIGS. 4 a and 4 b, calculating the separationfactor in terms of an “expansion factor”, k, by the simultaneous andequal expansion (k>1) or contraction (k<1) of both ellipses until theytouch, is neither too optimistic nor too pessimistic. This expansionfactor calculation can increase the allowable proximity between twoadjacent wells whilst satisfying the geometrical and probabilisticconstraints. Although iterative methods can be used starting from theelliptical conditions, there is no guarantee that such iterative schemesconverge towards the correct expansion factor solution.

SUMMARY

In embodiments of the invention, there is provided acomputer-implemented method and a system according to the appendedclaims.

According to an embodiment of the invention, there is provided acomputer-implemented method for determining the relative positions of awellbore and an object, the wellbore being represented by a firstellipse and the object being represented by a second ellipse, whereinthe first ellipse represents the positional uncertainty of the wellboreand the second ellipse represents the positional uncertainty of theobject, the method comprising the steps of:

receiving input data relating to a measured or estimated position of thewellbore and the object, the position of the wellbore having a first setof parameters defining the first ellipse, and the position of the objecthaving a second set of parameters defining the second ellipse;

calculating an expansion factor representing an amount by which one, orboth, of the first ellipse and the second ellipse can be expanded withrespect to one or both of respective first and second sets of ellipticalparameters so that the first and second ellipses osculate, whereincalculating the expansion factor involves determining and solving aquartic equation that is based on the geometry of the ellipses; and

determining, based on the calculated expansion factor, position dataindicative of the relative positions of the wellbore and the object.

Embodiments of the present invention utilize spatial relationshipsbetween two ellipses for determining the conditions governing osculationbetween the two ellipses (where osculation is the case in which theellipses touch), expressing the determination as an expansion scalefactor. Each expansion factor calculation involves using the smallestpositive root of the quartic equation. The explicit schemes of thepresent invention offer improvements in both calculation efficiency andreliability over known methods of calculating a separation factor andover iterative methods of calculating an expansion factor.

Typically, the wellbore is a first wellbore, and the object is a secondwellbore. Alternatively, the object may be a sub-surface hazard that isto be avoided when drilling the wellbore.

Methods are presented for the expansion of either one, or both ellipses.The computer-implemented methods can be used to increase the allowableproximity of two adjacent wellbores whilst satisfying the necessarygeometrical and probabilistic constraints. The calculation method isconsistent with existing industry wellbore uncertainty models. Since thedetermination of the osculating condition is exact the calculation isneither too optimistic nor too conservative.

Further features and advantages of the invention will become apparentfrom the following description of preferred embodiments of theinvention, given by way of example only, which is made with reference tothe accompanying drawings.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 shows a three-dimensional representation of a cone whichrepresents the positional uncertainty of a wellbore;

FIG. 2 shows a “centre vector method” for estimating the separationbetween two ellipses;

FIG. 3 shows a “pedal curve method” for estimating the separationbetween two ellipses;

FIG. 4 a shows the determination of an “expansion factor” by thesimultaneous and equal expansion (k>1) of two ellipses;

FIG. 4 b shows the determination of an “expansion factor” by thesimultaneous and equal contraction (k<1) of two ellipses;

FIG. 5 a shows the steps involved in a first method of calculating anexpansion factor;

FIGS. 5 b-5 e show an expansion of the ellipses carried out in the firstmethod of calculating the expansion factor;

FIG. 6 a shows the steps involved in a second method of calculating anexpansion factor;

FIG. 6 b shows an expansion of the ellipses carried out in the secondmethod of calculating the expansion factor;

FIGS. 7 a-9 b show dual and single sided expansion of variousconfigurations of ellipses;

FIG. 10 shows a wellbore positioning system according to the presentinvention;

FIG. 11 shows an example of a planned wellbore in simplified collisionavoidance plot;

FIG. 12 shows the steps taken in determining the relative position of awellbore according to the present invention;

FIG. 13 shows a schematic diagram of a wellbore being drilled into aformation.

DETAILED DESCRIPTION

In directional work an ellipse is generally defined by its centre (x₀,y₀), the lengths of its semi-major and semi-minor axes a and b and theorientation θ of the major axis direction a relative to some referencedirection.

Since the dimensions of each ellipse represents some confidence intervalthat the wellbore lies within its boundary, the equal expansion orcontraction of both ellipses until they touch (i.e. osculate) is ameasure of a potential collision between the wells. Since the point atwhich two ellipses touch is a function of both their sizes andorientations, conceptually, the available space can also be calculatedby expanding only one ellipse with the other one fixed. Therefore bothdual sided and single sided expansion can be applied to calculate arelevant expansion factor k.

Mathematically, it may be more convenient to represent the ellipse as aquadratic form, as shown by Eq. 5, incorporating the above ellipticalparameters and the expansion factor k within the quadratic form'scoefficients, which are shown by Eqs. 6 to 11. Details of both thetransform and inverse transform are given in Appendix A.

E(x,y,k)=Ax ²+2Bxy+Cy ²+2Dx+2Fy+H−a ² b ² k ²=0  (5)

Where

A=b ² cos² θ+a ² sin²θ  (6)

B=(b ² −a ²)sin θ cos θ  (7)

C=a ² cos² θ+b ² sin²θ  (8)

D=−y ₀ B−x ₀ A  (9)

F=−x ₀ B−y ₀ C  (10)

H=x ₀ ² A+2x ₀ y ₀ B+y ₀ ² C  (11)

The quadratic may be represented in matrix form, as shown by Eq. 12,where E is the 3×3 symmetric matrix. It is noted that the expansionfactor appears in only one of the matrix elements as its square k².

$\begin{matrix}{{E( {x,y,k} )} = {{\begin{bmatrix}x & y & 1\end{bmatrix}\begin{bmatrix}A & B & D \\B & C & F \\D & F & {H - {a^{2}b^{2}k^{2}}}\end{bmatrix}}\begin{bmatrix}x \\y \\1\end{bmatrix}}} & (12)\end{matrix}$

It is also possible to represent the ellipse in matrix form so that thesquare symmetric matrix is independent of the ellipse's origin, (Zheng,X., Palffy-Muhoray, P.: “Distance of Closest Approach of Two ArbitraryHard Ellipses in 2D”). A summary of the representation is included inAppendix A.

ZPM Expansion Factor

The following calculation, referred to hereinafter as the “ZPM” method,can be used to calculate the expansion factor using dual sidedexpansion, where each of two ellipses is expanded equally.

Referring to FIG. 5 a, in step S501 the elliptical parameters a₁, b₁,θ₁, x_(0,1), y_(0,1) of ellipse E₁ and a₂, b₂, θ₂, x_(0,2), y_(0,2) ofellipse E₂ are input into a wellbore positioning system, as describedbelow with respect to FIG. 10. In step S502, a determination is made asto whether the centres of the ellipses are separated by a distancegreater than δ_(min), as explained above in relation to Eq. 3. If theseparation is not greater than δ_(min), it is determined in step S503that the wellbores physically interfere and the calculation is stoppedin step S504. Alternatively, if the ellipse centres are separated bymore than δ_(min), the distance of closest approach δ_(cr) is calculatedin step S505, as explained below.

The distance of closest approach δ_(cr) of two arbitrary hard ellipsesin 2D can be determined using the method disclosed in Zheng, X.,Palffy-Muhoray, P.: “Distance of Closest Approach of Two Arbitrary HardEllipses in 2D”. Referring to FIG. 5 b, the ellipse E₂ is translatedtowards E₁ in the direction joining their centres until it reaches theposition E₂* when the ellipses touch externally. The orientations of thetwo ellipses are maintained throughout. The ellipse E₁ is thentransformed into a circle C₁ and the same mathematical transformationused to obtain the circle is applied to the ellipse E₂ (FIG. 5 c). Notethat the circle C₁ and the ellipse E₂** remain connected at a respectivetangent after the transformation. The closest approach between thecircle and the ellipse is then found analytically, recovering theclosest approach between the ellipses E₁ and E₂* by applying the inverseof the transformation used to obtain the circle (FIG. 5 d). Details ofthe relevant calculations are explained in the paper by Zheng andPalffy-Muhoray, together with the solution of the resulting quarticequation; the quartic equation is given in Appendix B.

The problem's symmetry is then used to determine the expansion factor(FIG. 5 e). It is noted that the translation of E₂ to E₂* followed by amagnification of magnitude k, of both E₁ and E₂* together (whilstmaintaining their relative position) about the centre of E₁ isequivalent to the magnification of magnitude k of each of E₁ and E₂about their respective centres, (Snapper, E., Troyer, R. J.: “MetricAffine Geometry”, 1971, Academic Press, London, 1, 36-55). Therefore,the dual sided expansion factor k can be calculated from the distance ofclosest approach using the scaling factor k=δ/δ_(cr), and the expansionfactor k is output (step S507).

Zheng and Palffy-Muhoray also describe a method for calculating thecontact point and provided computer code for both of the closestapproach and contact point calculations. Knowledge of the contact pointmay be used to verify the expansion factor results, checking for eachellipse that |E(x, y, k)|<ε, where ε is some acceptable tolerance.

YKC Expansion Factor

A further method, referred to hereinafter as the “YKC” method, can beused to calculate an expansion factor using dual sided expansion, whereeach of the two ellipses is expanded equally, or single sided expansion,where only one ellipse is expanded while the other remains fixed.However, for dual expansion the ZPM approach is preferred. Tests showthat it is more stable computationally, particularly for similarly sizedellipses with centres that are close together.

Referring to FIG. 6 a, in step S601 the elliptical parameters a₁, b₁,θ₁, x_(0,1), y_(0,1) of ellipse E₁ and a₂, b₂, θ₂, x_(0,2), y_(0,2) ofellipse E₂ are input into a wellbore positioning system, as describedbelow with respect to FIG. 10. In step S602, a determination is made asto whether dual sided or single sided expansion is preferred. Singlesided expansion may be preferred in some cases because of the greaterarea of space obtained about the expanded wellbore.

Single Sided Expansion

For single sided expansion (output “Y” at step S602), the symmetrypresent in the dual sided expansion is broken and a different approachmust be used. Referring to FIG. 6 a, in this case the size of the firstellipse E₁ is fixed and for a solution to exist, the centre (x_(0,2),y_(0,2)) of the second ellipse E₂ must lie outside its boundary (stepS603). Mathematically this requires the condition that E₁(x_(0,2),y_(0,2), 1)>0. In step S604, if the centre of E₂ does not lie outside ofE₁, the system determines that no solution is possible, and thecalculation is stopped at step S605.

A characteristic cubic polynomial P(λ)=det(λE ₁−E ₂)=0, which can beused to determine the separation conditions between two ellipses withoutexplicitly calculating the contact point, was derived in Choi, Y. K.:“Collision Detection for Ellipsoids and Other Quadrics”, PhD Thesis,University of Hong Kong, March 2008. Choi showed that if E₁ and E₂ aretwo ellipses with the characteristic polynomial P(λ) (where λ is amultiplier) then they are separated if and only if P(λ) has two distinctnegative roots and they touch each other externally if and only if P(λ)has a double negative root. The ellipses are overlapping if P(λ) has nonegative root.

For the purpose of calculating the expansion factor, the expansionfactor can be incorporated in the characteristic polynomial givingP(λ)=det[λE ₁(k₁)−E ₂(k₂)]=0. For a single sided expansion the firstellipse E₁ is fixed so set k₁=1 and k₂=k. The characteristic, cubicpolynomial becomes P(λ)=det[λE ₁−E ₂(k)]=0 (step S606). Using Choi'scondition, the cubic's discriminant vanishes when the ellipses touch,leaving a quartic equation in k², as shown by Eq.13. Taking the squareroot gives the expansion factor k.

γ₄ k ⁸+γ₃ k ⁶+γ₂ k ⁴+γ₁ k ²+γ₀=0  (13)

After lengthy, computer assisted simplification, using a softwareprogram such as Mathematica, the coefficients of the quartic equationcan be written as Eqs. 14 to 26. Further details are provided inAppendix B.

$\begin{matrix}{\mspace{79mu} {\gamma_{4} = {a_{2}^{4}{b_{2}^{4}( {z_{12}^{2} - {4a_{1}^{2}a_{2}^{2}b_{1}^{2}b_{2}^{2}}} )}}}} & (14) \\{\gamma_{3} = {2a_{2}^{2}{b_{2}^{2}\lbrack {{6a_{1}^{2}a_{2}^{2}b_{1}^{2}{b_{2}^{2}( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )}} + {9a_{1}^{2}a_{2}^{2}b_{1}^{2}b_{2}^{2}z_{12}a_{2}^{2}b_{2}^{2}r_{1}z_{12}} - {z_{12}^{2}( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )} - {2z_{12}^{3}}} \rbrack}}} & (15) \\{\gamma_{2} = {{{- 27}a_{1}^{4}b_{1}^{4}a_{2}^{4}b_{2}^{4}} + {a_{2}^{4}b_{2}^{4}r_{1}^{2}} - {6a_{1}^{2}b_{1}^{2}a_{2}^{2}b_{2}^{2} \times \lfloor {{2( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )^{2}} + {3a_{2}^{2}b_{2}^{2}r_{1}}} \rfloor} - {2a_{2}^{2}b_{2}^{2}{z_{12}( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )}( {{9a_{1}^{2}b_{1}^{2}} - {2r_{1}}} )} + {z_{12}^{2}( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )}^{2} + {12a_{2}^{2}b_{2}^{2}{r_{1}\lbrack {{( {{a_{2}^{2}b_{1}^{2}} + {a_{1}^{2}b_{2}^{2}}} ) \times {\cos^{2}( {\theta_{1} - \theta_{2}} )}} + {( {{a_{1}^{2}a_{2}^{2}} + {b_{1}^{2}b_{2}^{2}}} ){\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rbrack}^{2}}}} & (16) \\{\gamma_{1} = {2( {{\lbrack {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} \rbrack \{ {{{- a_{2}^{2}}b_{2}^{2}r_{1}^{2}} + {a_{1}^{2}{b_{1}^{2}\lbrack {{2( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )^{2}} + {9a_{2}^{2}b_{2}^{2}r_{1}}} \rbrack}}} \}} - {r_{1}z_{12}\lfloor {( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )^{2} + {6a_{2}^{2}b_{2}^{2}r_{1}}} \rfloor}} )}} & (17) \\{\mspace{79mu} {{\gamma_{0} = {r_{1}^{2}\lfloor {{4a_{2}^{2}b_{2}^{2}r_{1}} + ( {{a_{2}^{2}p_{1}} + {b_{2}^{2}q_{1}}} )^{2}} \rfloor}}\mspace{79mu} {Where}}} & (18) \\{\mspace{79mu} {\phi_{1} = {{\Delta \; y\; \cos \; \theta_{1}} - {\Delta \; x\; \sin \; \theta_{1}}}}} & (19) \\{\mspace{79mu} {\vartheta_{1} = {{\Delta \; x\; \cos \; \theta_{1}} + {\Delta \; y\; \sin \; \theta_{1}}}}} & (20) \\{\mspace{79mu} {\phi_{2} = {{\Delta \; y\; \cos \; \theta_{2}} - {\Delta \; x\; \sin \; \theta_{2}}}}} & (21) \\{\mspace{79mu} {\vartheta_{2} = {{\Delta \; x\; \cos \; \theta_{2}} + {\Delta \; y\; \sin \; \theta_{2}}}}} & (22) \\{\mspace{79mu} {p_{1} = {\phi_{2}^{2} - {a_{1}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} - {b_{1}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}}}} & (23) \\{\mspace{79mu} {q_{1} = {\vartheta_{2}^{2} - {a_{1}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} - {b_{1}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}}}} & (24) \\{\mspace{79mu} {r_{1} = {{b_{1}^{2}\vartheta_{1}^{2}} + {a_{1}^{2}( {\phi_{1}^{2} - b_{1}^{2}} )}}}} & (25) \\{z_{12} = {{a_{1}^{2}\lbrack {{a_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rbrack} + {b_{1}^{2}\lfloor {{a_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor}}} & (26)\end{matrix}$

By inspection, this calculation of the closest distance of a point(which may represent an object) to an ellipse is equivalent to thesingle sided expansion of a unit circle (which is a special case of anellipse) centred on the point against the ellipse, as shown in FIG. 6 b.This distance is equal to the expansion factor k (step S607). Oncecalculated, the expansion factor is output at step S608.

Dual Sided Expansion

For dual sided expansion, in step S609 a determination is made as towhether the centres of the ellipses are separated by a distance greaterthan δ_(min), as explained above in relation to Eq. 3. If the separationis not greater than δ_(min), it is determined in step S610 that thewellbores physically interfere and the calculation is stopped in stepS611. Alternatively, if the ellipse centres are separated by more thanδ_(min), dual sided expansion can proceed; both ellipses are expandedequally so set k₁=k₂=k. The characteristic polynomial becomesP(λ)=det[λE ₁(k)−E ₂(k)]=0 (step S612). This results in another quarticequation in k²; details of the coefficients are provided in Appendix B.When solved (step S613), the smallest positive root of the equationgives the expansion factor k, at which point the calculation stops (stepS614).

EXAMPLES

Some examples of elliptical configurations are shown in FIGS. 7 a to 9b. The configurations of two ellipses on FIGS. 7 a, 8 a and 9 acorrespond to the configurations in FIGS. 7 b, 8 b and 9 b,respectively. In FIGS. 7 a, 8 a and 9 a, the dashed ellipses representthe expanded, osculating ellipses when a dual expansion method is used.FIGS. 7 b, 8 b and 9 b the dashed ellipse represents the expansion ofone of the ellipses in a single sided expansion.

In using the ZPM dual sided expansion method and the YKC single sidedexpansion method, the expansion factors for these configurations arecalculated as follows:

FIG. 7 a—the dual sided expansion for the ellipses E₁(3, 2, 0°, 0, 0)and E₂(4, 2, 90°, 8, 0) gives the expansion factor k=1.6;

FIG. 7 b—the single sided expansion for the ellipses E₁(3, 2, 0°, 0, 0)and E₂(4, 2, 90°, 8, 0) gives the expansion factor k=2.5;

FIG. 8 a—the dual sided expansion for the ellipses E₁(3, 2, 30°, −3, −2)and E₂(2, 1, 135°, 1, 1) gives the expansion factor k=1.25568;

FIG. 8 b—the single sided expansion for the ellipses E₁(3, 2, 30°, −3,−2) and E₂(2, 1, 135°, 1, 1) gives the expansion factor k=2.01033;

FIG. 9 a—the dual sided expansion for the ellipses E₁(7, 1, 135°, −2, 0)and E₂(5, 1, 150°, 1, 0) gives the expansion factor k=0.814767; and

FIG. 9 b—the single sided expansion for the ellipses E₁(7, 1, 135°, −2,0) and E₂(5, 1, 150°, 1, 0) gives the expansion factor k=0.695637.

Referring to the configurations of FIGS. 7 a, 8 a and 9 a, Table 1 showsa comparison of the CVM and PCM separation factors (k_(CVM) and k_(PCM),respectively) with the dual sided expansion factor (k_(ZPM)) for thethree elliptical configurations.

Factor Ellipses k_(CVM) k_(PCM) k_(ZPM) E₁ E₂ [−] [−] [−] 3, 2, 0°, 0, 04, 2, 90°, 8, 0 1.60000 1.60000 1.60000 3, 2, 30°, −3, −2 2, 1, 135°, 1,1 1.25593 1.24452 1.25568 7, 1, 135°, −2, 0 5, 1, 150°, 1, 0 0.911900.32055 0.81477

From the separation and expansion factors of Table 1, it can be seenthat the factors are calculated to be the same value (1.6) for theconfiguration of FIG. 7 a; calculations agree only in a special casewhere the major or minor axes of the ellipses are collinear. Thedifferences in the calculated factors in any particular case may be muchmore pronounced as eccentricities increase.

System

As shown in FIG. 13, wellbore drilling systems generally comprisedrilling equipment 4 arranged to drill a wellbore 5 into the one or morehydrocarbon-bearing reservoirs in a formation 6. The drilling systemtypically comprises a controller 7 arranged to control the drillingequipment. An existing wellbore 8 is also shown.

In order to determine optimum settings of the various components of thewellbore drilling system, the wellbore positioning system 100 comprisessuitable computer-implemented models, software tools and hardware, asshown in FIG. 10. A reservoir model 121 may be employed. As known in theart, a reservoir model is a conceptual 3-dimensional construction of areservoir that is constructed from incomplete data with much of theinter-well space estimated from data obtained from nearby wells or fromseismic data. In conjunction with this, a trajectory model 123, that is,a computer model that constructs 2D and/or 3D representations of thegeographical locations and/or trajectories of wellbores may be employed.The trajectory model may comprise or make use of a collision avoidanceplot, also known as a travelling cylinder or normal plane diagram. Anexpansion factor calculation tool 111, as described further below, cancalculate the expansion factor as explained above. Using the reservoirmodel 121, the trajectory model 123 can use information such as thevolume and shape of the reservoir 3 (including the arrangement ofoverlying rock formations and the locations of any faults or fracturesin the rock formations and sub-surface hazards), the porosity of theoil-bearing rock formations, the location of existing production well(s)and injection well(s), in combination with the results of the expansionfactor calculation tool 111, to provide an indication as to the possibletrajectory of a planned wellbore.

In one arrangement, referring to FIG. 10, the expansion factorcalculation tool 111 and optionally the reservoir model 121, thetrajectory model 123 and an optimisation tool 125 are executed by thewellbore positioning system 100. The wellbore positioning system 100,which is for example a control system on a platform, can compriseconventional operating system and storage components such as a systembus connecting a central processing unit (CPU) 105, a hard disk 103, arandom access memory (RAM) 101, and I/O and network adaptors 107facilitating connection to user input/output devices and interconnectionwith other devices on a network N1. The Random Access Memory (RAM) 101contains operating system software 131 which controls, in a knownmanner, low-level operation of the wellbore positioning system 100. Theserver RAM 101 contains the software tools and models 111, 121, 123 and125 during execution thereof. Each item of software is configurable withmeasurement and/or predetermined data stored in a database or otherstorage component which is operatively coupled or connected to thewellbore positioning system 100; in the system of FIG. 2, storagecomponent DB1 stores all such data relating to the expansion factorcalculation tool 111 and is accessible thereby, while storage componentDB2 stores all other data for use by the other components of the system100.

Input data received by receiving means of the system 100 comprise theelliptical parameter values and are based on a measured position of anexisting wellbore or an estimated (i.e. modelled or simulated) positionof a planned wellbore. Such estimated input data can be modelled orestimated upon planning a wellbore, for example upon an initialassessment or appraisal of a reservoir when developing a new field.Alternatively, in the case where the position of a planned wellbore isbeing determined in order to avoid an object other than anotherwellbore, such as a sub-surface hazard, the input data includesmeasurement data relating to the position of the object.

The measurement data may comprise specific measured values as directlymeasured by suitably positioned measurement equipment such as surveyinstruments 12, or may comprise values derived from a number of separatepositional measurements. Therefore, the raw measured data may, ifnecessary or preferred, be manipulated by appropriate software andexecuted by the CPU 105 of the system 100, in order to generatemeasurement or estimated position data that are suitable for inputtinginto the expansion factor calculation tool 111. Such manipulation maycomprise using the reservoir and/or trajectory models to determine theparameter values of the two ellipses.

The expansion factor calculation tool 111 may comprise a softwareprogram such as Mathematica. This program can be used in a number ofways during the calculation of the expansion factor. Firstly by makinguse of its symbolic manipulation, the substitutions, for example, for A,B, C, D, F, G (which is equivalent to H−a²b²k²—see Appendix A), H can bemade. The determinants can then be expanded and the equations simplifiedusing this program. Additionally, Mathematica™ is preferably employed toprogram the resulting quartic coefficients and solve the quarticequation. Alternatively, the expressions can be programmed in, forexample, Visual Basic™ within an EXCEL™ spreadsheet.

An optimisation tool 125 may be provided to assist in the planning anddrilling of wellbores. The optimisation tool may be used in conjunctionwith the trajectory model 123 to compute an optimal position for thewellbore in 2D or an optimal trajectory in 3D, based on input dataincluding the calculated expansion factor and the measured or estimatedinput data that relates to the position of one or more existingwellbores or objects. In the case where a number of positions ortrajectories are possible, the optimisation tool 125 may be programmedwith rules that take into account additional data representing, forexample, threshold values representing practical limits to the degree ofcurvature of the wellbore trajectory. In this way, the optimisation tool125 can determine an optimum alignment of the trajectory, as explainedfurther below with reference to FIG. 11.

FIG. 11 shows a simplified collision avoidance plot which may beproduced by the trajectory model 123 upon calculation of the expansionfactor; the x and y axes represent length in metres. In FIG. 11 thedashed ellipse represents the tolerable errors, including an acceptableoperational margin, for a planned wellbore at some point in space. Thesolid ellipses represent the tolerable errors surrounding threeadjacent, drilled wellbores. By inspection, at this position in thewellbore the planned wellbore is heavily constrained and its position atthis point cannot be moved within the collision avoidance plot withoutinfringing the space in which the other wellbores may lie. The centrevector method is generally excluded in such a scenario as it is overlyoptimistic. Using the pedal curve method, a well planner would concludethat the planned wellbore could not be threaded through this point. Awell planner using the expansion factor calculation method, whichhonours the geometry, would conclude that the well could, with care,pass through this point. This is confirmed by the common sense approachthat, visually, the dashed ellipse fits comfortably within the availablespace.

The use of the expansion factor in the wellbore positioning method andsystem of the invention is advantageous in the planning and drilling ofwellbores, as it provides more space in which to plan and optimise thetrajectories of wellbores. However, if a planner concluded that it wasnot possible to drill through the gap of FIG. 11, then the wellborewould have to be planned around the existing wellbores. Such activitiesadd to the tortuosity of the wellbore's trajectory, which increasestorque and drag forces, and/or may be difficult to achieve with theavailable tools. In some cases the detour may not be possible. Insubsurface terms, the detour may make it difficult to achieve optimumalignment to a target. If so, oil and gas reserves and production may beadversely affected.

The wellbore positioning system 100 is preferably operatively connectedto a controller 133 of the wellbore drilling system, for example via thenetwork N1. The controller 133 of the wellbore drilling system isautomatically configured with the one or more operating modes determinedby the system 100, the controller 133 being arranged to apply the one ormore operating modes.

Method

Referring to FIG. 12, the steps involved in a first embodiment of acomputer-implemented method for determining one or more operating modesfor the wellbore drilling system are shown.

In step S1201, the input data is received by the wellbore positioningsystem 100.

At step S1202, the input data are input into the expansion factorsoftware tool 111, the calculations of which are described above inrelation to FIGS. 5 a-5 e, 6 a and 6 b. The expansion factor calculationtool is then run in step S1203, and generates, at step S1204, positiondata indicative of a relative position or proximity of the plannedwellbore to the existing or simulated wellbore or object. This data maybe output in various forms, for example, as coordinates of a 2D or 3Dsimulation of a reservoir, or as a collision avoidance plot.

At step S1205, the generated position data are used to determine one ormore operating modes of the wellbore drilling system. The operating modecan represent an instruction or suggested setting for the drillingsystem, which can subsequently be applied to the drilling system. Thedetermination can include the step of comparing, in accordance with apredetermined set of rules (which can be set using a collision avoidanceplot implemented by the trajectory model 123), the calculated positiondata to predetermined known or threshold position data that isaccessible from the database DB2. For example, the determination may bebased on a known position of an existing wellbore or a sub-surfacehazard.

Software executed by the CPU 105 of the system 100 determines, on thebasis of the determined position data, the one or more operating modesof the wellbore drilling system. The expansion factor calculation tool111, the reservoir model 121 and/or the trajectory model 123 may beconfigured to determine the operating mode(s) upon generation of theposition data, or a separate software component may be provided.Additional technical and physical constraints determined by thereservoir model 121 or the trajectory model 123 may be taken intoaccount in order to determine the operating mode, and can be stored andaccessed from the databases DB1 and DB2 as necessary.

For example, the operating mode can comprise an instruction to go aheadwith the drilling of a planned wellbore or not, this determination beingbased on a determination by the trajectory model 123 that the trajectoryof the planned wellbore under consideration is drillable. Alternativelyor additionally, the operating mode can comprise one or more specificconfiguration settings for the wellbore drilling system, such as adrilling speed or trajectory.

The software component used to determine the operating mode isconfigured to use a predetermined set of rules in conjunction with inputdata such as the calculated expansion factor, in order to determine theoperating mode. These rules are stored in and accessible from thedatabase DB1 and DB2 as necessary.

The computer-implemented method can further include an optional step,S1206, of applying or inputting the determined operating mode into acontroller of the wellbore drilling system.

The above embodiments are to be understood as illustrative examples ofthe invention. It is to be understood that any feature described inrelation to any one embodiment may be used alone, or in combination withother features described, and may also be used in combination with oneor more features of any other of the embodiments, or any combination ofany other of the embodiments. Furthermore, equivalents and modificationsnot described above may also be employed without departing from thescope of the invention, which is defined in the accompanying claims.

APPENDIX A Ellipse Representations

The derivations using the YKC conditions depend on the ability totranslate freely between the ellipse representations. The first andsecond quadratic forms are mathematically equivalent.

First Quadratic Form

The ellipse E₀(x, y)=b²x²+a²y²−a²b²=0 with semi-major axis a andsemi-minor axis b, aligned with the x and y axes and centred on theorigin can be represented as a quadratic form, Eq. A-1.

$\begin{matrix}{{E_{0}( {x,y} )} = {{\begin{bmatrix}x & y & 1\end{bmatrix}\begin{bmatrix}b^{2} & 0 & 0 \\0 & a^{2} & 0 \\0 & 0 & {{- a^{2}}b^{2}}\end{bmatrix}}\begin{bmatrix}x \\y \\1\end{bmatrix}}} & ( {A\text{-}1} )\end{matrix}$

Writing the ellipse E₀(x, y)=x E ₀ x ^(T) gives

$\begin{matrix}{{\underset{\_}{E}}_{0} = \begin{bmatrix}b^{2} & 0 & 0 \\0 & a^{2} & 0 \\0 & 0 & {{- a^{2}}b^{2}}\end{bmatrix}} & ( {A\text{-}2} )\end{matrix}$

The matrix T translates a point on the ellipse by an amount x₀ in the xdirection and y₀ in the y direction. The rotation matrix R rotates apoint by an amount θ clockwise about the origin. The scaling matrix Sscales a point by a factor k relative to the origin.

$\begin{matrix}{\underset{\_}{T} = \begin{bmatrix}1 & 0 & 0 \\0 & 1 & 0 \\{- x_{0}} & {- y_{0}} & 1\end{bmatrix}} & ( {A\text{-}3} ) \\{\underset{\_}{R} = \begin{bmatrix}{\cos \; \theta} & {{- \sin}\; \theta} & 0 \\{\sin \; \theta} & {\cos \; \theta} & 0 \\0 & 0 & 1\end{bmatrix}} & ( {A\text{-}4} ) \\{\underset{\_}{S} = \begin{bmatrix}k^{- 1} & 0 & 0 \\0 & k^{- 1} & 0 \\0 & 0 & 1\end{bmatrix}} & ( {A\text{-}5} )\end{matrix}$

Note that the translation, rotation and scaling operations do notcommute and therefore the order in which these operations are performedis important. Combining these transformations as shown in Eq. A-6 givesthe ellipse E(x, y, k) in the body of the paper, Eqs. 5 and 12.

E=x T R S E ₀ S ^(T) R ^(T) T ^(T) x ^(T)  (A-6)

Inverse Transform

The coordinates of the ellipse's centre, semi-major and semi-minor axesand orientation may be recovered from the first quadratic form using theinverse transform, Eqs. A-7 to A-11. Although the inverse transform isnot used in either the separation or expansion factor calculations itprovides an effective means of testing the correctness of the transform.Note that the constant G is equivalent to H−a²b²k².

$\begin{matrix}{x_{0} = \frac{{CD} - {BF}}{B^{2} - {AC}}} & ( {A\text{-}7} ) \\{y_{0} = \frac{{AF} - {BD}}{B^{2} - {AC}}} & ( {A\text{-}8} ) \\{a = \sqrt{\frac{2( {{AF}^{2} + {CD}^{2} + {GB}^{2} - {2{BDF}} - {ACG}} )}{( {B^{2} - {AC}} )\lbrack {\sqrt{( {A - C} )^{2} + {4B^{2}}} - ( {A + C} )} \rbrack}}} & ( {A\text{-}9} ) \\{b = \sqrt{\frac{2( {{AF}^{2} + {CD}^{2} + {GB}^{2} - {2{BDF}} - {ACG}} )}{( {B^{2} - {AC}} )\lbrack {{- \sqrt{( {A - C} )^{2} + {4B^{2}}}} - ( {A + C} )} \rbrack}}} & ( {A\text{-}10} ) \\{\theta = \{ \begin{matrix}0 & {{B = 0},{A < C}} \\{\pi/2} & {{B = 0},{A > C}} \\{\frac{1}{2}{\tan^{- 1}( \frac{2B}{A - C} )}} & {{B \neq 0},{A < C}} \\{\frac{\pi}{2} + {\frac{1}{2}{\tan^{- 1}( \frac{2B}{A - C} )}}} & {{B \neq 0},{A > C}}\end{matrix} } & ( {A\text{-}11} )\end{matrix}$

Second Quadratic Form

The ellipse can also be represented so the symmetric matrix isindependent of the ellipse's origin, Eq. A-12 and A-13, (Zheng andPalffy-Muhoray, 2010). Here I is the identity matrix and the vectorθ=[sin θ, cos θ]. In Zheng and Palffy-Muhoray's paper these authorsassume k=1 throughout.

$\begin{matrix}{{{{\lbrack {{x - x_{0}},{y - y_{0}}} \rbrack \begin{bmatrix}m_{11} & m_{12} \\m_{21} & m_{22}\end{bmatrix}}\begin{bmatrix}{x - x_{0}} \\{y - y_{0}}\end{bmatrix}} = k^{2}}{Where}} & ( {A\text{-}12} ) \\{\begin{bmatrix}m_{11} & m_{12} \\m_{21} & m_{22}\end{bmatrix} = {\frac{1}{b^{2}}\lbrack {\underset{\_}{I} + {( {\frac{b^{2}}{a^{2}} - 1} ){\underset{\_}{\theta}}^{T}\underset{\_}{\theta}}} \rbrack}} & ( {A\text{-}13} )\end{matrix}$

APPENDIX B Expansion Factors YKC Expansion Factor (Single Sided)

For a single sided expansion the characteristic polynomial becomesP(λ)=det[λE ₁−E ₂(k)]=0, Eq. B-1. For conciseness substitute χ=k².Expanding the determinant gives a cubic polynomial which coefficientsare functions of the coefficients of the quadratic forms and the squareof the expansion factor, Eq. B-2.

$\begin{matrix}{{\det \{ {{\lambda \begin{bmatrix}A_{1} & B_{1\;} & D_{1} \\B_{1} & C_{1} & F_{1} \\D_{1} & F_{1} & {H_{1} - {a_{1}^{2}b_{1}^{2}}}\end{bmatrix}} - \begin{bmatrix}A_{2} & B_{2} & D_{2\;} \\B_{2} & C_{2\;} & F_{2} \\D_{2} & F_{2} & {H_{2} - {a_{2}^{2}b_{2}^{2}\chi}}\end{bmatrix}} \}} = 0} & ( {B\text{-}1} ) \\{\mspace{79mu} {{{{w_{3}(\chi)}\lambda^{3}} + {{w_{2}(\chi)}\lambda^{2}} + {{w_{1}(\chi)}\lambda} + {w_{0}(\chi)}} = 0}} & ( {B\text{-}2} )\end{matrix}$

Choi, 2008 showed that the cubic discriminant equals zero when theellipses touch externally, Eq. B-3.

w ₂ ² w ₁ ²−4w ₂ ³ w ₀+18w ₃ w ₂ w ₁ w ₀−4w ₃ w ₁ ³−27w ₃ ² w ₀²=0  (B-3)

By inspection, the coefficients w_(j) are of the form w_(j)=u_(j)+v_(j)χ

Making this substitution, Eq. B-3 then becomes a quartic equation in χ,Eq. B-4.

$\begin{matrix}{{{( {{v_{1}^{2}v_{2}^{2}} - {4v_{1}^{3}v_{3}}} ){\chi^{4}( {{2u_{2}v_{1}^{2}v_{2}} + {2u_{1}v_{1}v_{2}^{2}} - {4u_{0}v_{2}^{3}} - {12u_{1}v_{1}^{2}v_{3}} + {18u_{0}v_{1}v_{2\;}v_{3}}} )}\chi^{3}} + {( {{u_{2}^{2}v_{1}^{2}} + {4u_{1}u_{2}v_{1}v_{2}} + {u_{1}^{2}v_{2}^{2}} - {12u_{0}u_{2}v_{2}^{2}} - {12u_{1}^{2}v_{1}v_{3}} + {18u_{0}u_{2}v_{1}v_{3}} + {18u_{0}u_{1}v_{2}v_{3}} - {27u_{0}^{2}v_{3}^{2}}} )\chi^{2}} + {( {{2u_{1}u_{2}^{2}v_{1}} + {2u_{1}^{2}u_{2}v_{2}} - {12u_{0}u_{2}^{2}v_{2}} - {4u_{1}^{3}v_{3}} + {18u_{0}u_{1}v_{3}}} )\chi} + {u_{1}^{2}u_{2}^{2}} - {4u_{0}u_{2}^{3}}} = 0} & ( {B\text{-}4} )\end{matrix}$

Note that both u₃=0 and v₀=0 and

$\begin{matrix}{\mspace{79mu} {u_{0} = {{- a_{1}^{4}}b_{1}^{4}}}} & ( {B\text{-}5} ) \\{u_{1} = {a_{1}^{2}b_{1}^{2}\{ {{{- a_{2}^{2}}b_{1}^{2}{\cos^{2}( {\theta_{2} - \theta_{1}} )}} - {a_{1}^{2}\lfloor {{a_{2}^{2}{\sin^{2}( {\theta_{2} - \theta_{1}} )}} + {b_{2}^{2}{\cos^{2}( {\theta_{2} - \theta_{1}} )}}} \rfloor} + {b_{2}^{2}\lfloor {{{- b_{1}^{2}}{\sin^{2}( {\theta_{2} - \theta_{1}} )}} + ( {{\Delta \; x_{0}\cos \; \theta_{2}} + {\Delta \; y_{0}\sin \; \theta_{2}}} )^{2}} \rbrack} + {a_{2}^{2}( {{\Delta \; y_{0}\cos \; \theta_{2}} - {\Delta \; x_{0}\sin \; \theta_{2}}} )}^{2}} \}}} & ( {B\text{-}6} ) \\{u_{2} = {a_{2}^{2}{b_{2}^{2}\lbrack {{{- a_{1}^{2}}b_{1}^{2}} + {b_{1}^{2}( {{\Delta \; x_{0}\cos \; \theta_{1}} + {\Delta \; y_{0}\sin \; \theta_{1}}} )}^{2} + {a_{1}^{2}( {{\Delta \; y_{0}\cos \; \theta_{1}} - {\Delta \; x_{0}\sin \; \theta_{1}}} )}^{2}} \rbrack}}} & ( {B\text{-}7} ) \\{\mspace{79mu} {v_{1} = {{- a_{1}^{2}}b_{1}^{2}a_{2}^{2}b_{2}^{2}}}} & ( {B\text{-}8} ) \\{v_{2} = {{- a_{2}^{2}}b_{2}^{2}\{ {{a_{1}^{2}\lfloor {{a_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor} + {b_{1}^{2}\lfloor {{a_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor}} \}}} & ( {B\text{-}9} ) \\{\mspace{79mu} {v_{3} = {{- a_{2}^{2}}b_{2}^{2}}}} & ( {B\text{-}10} )\end{matrix}$

Then substitute Eqs. B-5 to B-10 into Eq. B-4. Simplification gives thecoefficients γ_(i) of the quartic equation, Eq. 13 given in the specificdescription.

YKC Expansion Factor (Dual Sided)

The derivation of the dual sided expansion with the characteristicpolynomial P(λ)=deq[λE ₁(k)−E ₂(k)]=0 proceeds in the same way. Againfor conciseness substitute χ=k².

$\begin{matrix}{{\det \{ {{\lambda \begin{bmatrix}A_{1} & B_{1\;} & D_{1} \\B_{1} & C_{1} & F_{1} \\D_{1} & F_{1} & {H_{1} - {a_{1}^{2}b_{1}^{2}\chi}}\end{bmatrix}} - \begin{bmatrix}A_{2} & B_{2} & D_{2\;} \\B_{2} & C_{2\;} & F_{2} \\D_{2} & F_{2} & {H_{2} - {a_{2}^{2}b_{2}^{2}\chi}}\end{bmatrix}} \}} = 0} & ( {B\text{-}11} )\end{matrix}$

The coefficients γ_(i) of the quartic equation may be calculated as Eq.B-12 to B-24.

$\begin{matrix}{\mspace{79mu} {\phi_{1} = {{\Delta \; y\; \cos \; \theta_{1}} - {\Delta \; x\; \sin \; \theta_{1}}}}} & ( {B\text{-}12} ) \\{\mspace{79mu} {\vartheta_{1} = {{\Delta \; x\; \cos \; \theta_{1}} + {\Delta \; y\; \sin \; \theta_{1}}}}} & ( {B\text{-}13} ) \\{\mspace{79mu} {\phi_{2} = {{\Delta \; y\; \cos \; \theta_{2}} - {\Delta \; x\; \sin \; \theta_{2}}}}} & ( {B\text{-}14} ) \\{\mspace{79mu} {\vartheta_{2} = {{\Delta \; x\; \cos \; \theta_{2}} + {\Delta \; y\; \sin \; \theta_{2}}}}} & ( {B\text{-}15} ) \\{\mspace{79mu} {r_{1} = {{\phi_{1}^{2}a_{1}^{2}} + {\vartheta_{1}^{2}b_{1}^{2}}}}} & ( {B\text{-}16} ) \\{\mspace{79mu} {r_{2} = {{\phi_{2}^{2}a_{2}^{2}} + {\vartheta_{2}^{2}b_{2}^{2}}}}} & ( {B\text{-}17} ) \\{p_{1} = {{a_{1}^{2}\lfloor {{a_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + b_{1}^{2} + {b_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor} + {b_{1}^{2}\lfloor {{a_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor}}} & ( {B\text{-}18} ) \\{p_{2} = {{a_{2}^{2}\lfloor {{a_{1}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{1}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor} + {b_{2}^{2}\lfloor {{a_{1}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} + a_{2}^{2} + {b_{1}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rfloor}}} & ( {B\text{-}19} ) \\{\gamma_{4} = {\frac{1}{2}\begin{Bmatrix}{{b_{1}^{2}b_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + {a_{2}^{2}\lbrack {{b_{1}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} - b_{2}^{2}} \rbrack} +} \\{a_{1}^{2}\lbrack {{a_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} - b_{1}^{2} + {b_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rbrack}\end{Bmatrix}^{2} \times ( {{2{a_{1}^{4}\lbrack {{a_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \rbrack}^{2}} + {2{b_{1}^{4}\lbrack {{a_{2}^{2}{\cos^{2}( {\theta_{1} - \theta_{2}} )}} + {b_{2}^{2}{\sin^{2}( {\theta_{1} - \theta_{2}} )}}} \rbrack}^{2}} + {a_{1}^{2}b_{1}^{2}\{ {{a_{2}^{2}{b_{2}^{2}\lbrack {{- 5} + {\cos \; 4( {\theta_{1} - \theta_{2}} )}} \rbrack}} + {4( {a_{2}^{4} + b_{2}^{4}} ){\sin^{2}( {\theta_{1} - \theta_{2}} )}{\cos^{2}( {\theta_{1} - \theta_{2}} )}}} \}}} )}} & ( {B\text{-}20} ) \\{\gamma_{3} = {2( {{{- 9}a_{1}^{2}b_{1}^{2}a_{2}^{2}b_{2}^{2}r_{1}p_{2}} + {6a_{1}^{2}b_{1}^{2}r_{2}p_{2}^{2}} - {9a_{1}^{2}b_{1}^{2}a_{2}^{2}b_{2}^{2}r_{2}p_{1}} - {r_{1}p_{2}^{2}p_{1}} + {6a_{2}^{2}b_{2}^{2}r_{1}p_{1}^{2}} - {r_{2}p_{2}p_{1}^{2}}} )}} & ( {B\text{-}21} ) \\{\gamma_{2} = {{18a_{1}^{2}b_{1}^{2}a_{2}^{2}b_{2}^{2}r_{1}r_{2}} - {12a_{1}^{2}b_{1}^{2}r_{2}^{2}p_{2}} + {r_{1}^{2}p_{2}^{2}} - {12a_{2}^{2}b_{2}^{2}r_{1}^{2}p_{1}} + {4r_{1}r_{2}p_{1}p_{2}} + {r_{2}^{2}p_{1}^{2}}}} & ( {B\text{-}22} ) \\{\mspace{79mu} {\gamma_{1} = {2( {{2a_{2}^{2}b_{2}^{2}r_{1}^{3}} + {2a_{1}^{2}b_{1}^{2}r_{2}^{3}} - {r_{1}^{2}r_{2}p_{2}} - {r_{2}^{2}r_{1}p_{1}}} )}}} & ( {B\text{-}23} ) \\{\mspace{79mu} {\gamma_{0} = {r_{1}^{2}r_{2}^{2}}}} & ( {B\text{-}24} )\end{matrix}$

ZPM Expansion Factor

Zheng and Palffy-Muhoray approach leads to a quartic equation in thevariable Q, Eq. (B-25).

$\begin{matrix}{{\tan^{2}{\varphi ( {\zeta + 1 - Q^{2}} )}( {\frac{Q}{b_{2}^{\prime}} + 1} )^{2}} = {( {Q^{2} - 1} )( {\frac{Q}{b_{2}^{\prime}} + 1 + \zeta} )^{2}}} & ( {B\text{-}25} )\end{matrix}$

This can be written in the standard form, Eq. (B-26)

ψ₄ Q ⁴+ψ₃ Q ³+ψ₂ Q ²+ψ₁ Q+ψ ₀=0  (B-26)

Where the quartic coefficients ψ_(i) are

$\begin{matrix}{\psi_{4} = {{- \frac{1}{b_{2}^{\prime 2}}}( {1 + {\tan^{2}\varphi}} )}} & ( {B\text{-}27} ) \\{\psi_{3} = {{- \frac{2}{b_{2}^{\prime 2}}}( {1 + \zeta + {\tan^{2}\varphi}} )}} & ( {B\text{-}28} ) \\{\psi_{2} = {{{- \tan^{2}}\varphi} - ( {1 + \zeta} )^{2} + {\frac{1}{b_{2}^{\prime 2}}\lbrack {1 + {( {1 + \zeta} )\tan^{2}\varphi}} \rbrack}}} & ( {B\text{-}29} ) \\{\psi_{1} = {\frac{2}{b_{2}^{\prime}}( {1 + {\tan^{2}\varphi}} )( {1 + \zeta} )}} & ( {B\text{-}30} ) \\{\psi_{0} = {( {1 + \zeta + {\tan^{2}\varphi}} )( {1 + \zeta} )}} & ( {B\text{-}31} )\end{matrix}$

To avoid a clash of symbols, note that the nomenclature used herediffers from that used by Zheng and Palffy-Muhoray. The variables ζ andφ used here are defined in their paper.

APPENDIX C Nomenclature

a=Ellipse semi-major axis length, L, ft

a=Unit vector in the major axis direction

b=Ellipse semi-minor axis length, L, ft

d=Diameter, L, ft

k=Expansion scale factor, dimensionless

m=Element of a transformation matrix

p=Substituted variable

r=Substituted variable

s=Characteristic length for a separation factor, L, ft

u=Substituted variable

v=Substituted variable

w=Coefficients of the YKC cubic equation

x=Ordinate in the normal plane, L, ft

y=Ordinate in the normal plane, L, ft

z=Substituted variable

A=First ellipse quadratic form coefficient

B=Second ellipse quadratic form coefficient

C=Third ellipse quadratic form coefficient

D=Fourth ellipse quadratic form coefficient

E=Ellipse

F=Fifth ellipse quadratic form coefficient

G=Sixth ellipse quadratic form coefficient

H=Modified sixth ellipse quadratic form coefficient

E=Ellipse matrix representation

R=Rotation matrix

S=Scaling matrix

T=Translation matrix

I=Unit matrix

P=Polynomial

Q=Independent variable of the ZPM quartic equation

Greek Symbols

δ=Centre to centre distance between ellipses, L, ft

γ=Coefficients of the YKC quartic equation

φ=Substituted or temporary variable

∂=Substituted or temporary variable

ψ=Coefficients of the ZPM quartic equation

χ=Square of the expansion factor, dimensionless

θ=Ellipse orientation angle to major axis, radians

φ=First variable defined by ZPM=

ζ=Second variable defined by ZPM

λ=Multiplier

Δ=A difference in a parameter

Subscripts and Superscripts

0=Condition at an origin

′=Transformed condition

1,2,3,4=First, second etc.

c=Casing

cr=Critical, or closest approach

h=Hole

i=Index

j=Index

s=Separation factor

1. A computer-implemented method for determining the relative positionsof a wellbore and an object, the wellbore being represented by a firstellipse and the object being represented by a second ellipse, whereinthe first ellipse represents a positional uncertainty of the wellboreand the second ellipse represents a positional uncertainty of theobject, the method comprising the steps of: receiving input datarelating to a measured or estimated position of the wellbore and theobject, the position of the wellbore having a first set of parametersdefining the first ellipse, and the position of the object having asecond set of parameters defining the second ellipse; calculating anexpansion factor representing an amount by which one, or both, of thefirst ellipse and the second ellipse can be expanded with respect to oneor both of respective first and second sets of elliptical parameters sothat the first and second ellipses osculate, wherein calculating theexpansion factor comprises determining and solving a quartic equationthat is based on the geometry of the first and second ellipses; anddetermining, based on the calculated expansion factor, position dataindicative of the relative positions of the wellbore and the object. 2.The method of claim 1, wherein the first and second ellipses areexpanded equally, and wherein the calculation of the expansion factorfurther comprises: solving the quartic equation to determine a distancebetween a centre of the first ellipse and a centre of the second ellipsewhen the second ellipse is translated towards the first ellipse along aline joining the centres of the first and second ellipses so that thefirst and second ellipses osculate; and calculating the expansion factorbased on the determined distance and a scale factor.
 3. The method ofclaim 1, wherein either the first and second ellipses are expandedequally, or one of the first and second ellipses is expanded, so thatthe first and second ellipses osculate, and wherein the calculation ofthe expansion factor further comprises: applying the first and secondsets of elliptical parameters to a polynomial equation, the solution ofwhich represents a separation condition of the ellipses; determining thequartic equation from the polynomial equation; and solving the quarticequation to calculate the expansion factor.
 4. The method of claim 3,wherein the ellipses osculate when the polynomial equation has a doubleroot.
 5. The method of claim 1, wherein the wellbore is a first plannedor drilled wellbore, and the object is a second planned or drilledwellbore.
 6. The method of claim 1, wherein the object is a sub-surfacehazard.
 7. The method of claim 1, wherein the first and second sets ofelliptical parameters are derived from a measured or estimated positionof the wellbore and the object.
 8. The method of claim 1, furthercomprising the step of determining a trajectory of the wellbore in athree-dimensional simulation.
 9. The method of claim 1, furthercomprising the step of optimising the position of the wellbore relativeto the object.
 10. A computer-implemented method for determining one ormore operating modes of a wellbore drilling system, the wellboredrilling system being arranged to drill a wellbore in a rock formation,the method comprising the steps of: receiving position data determinedaccording to the method of claim 1; inputting said position data into awellbore trajectory model; operating the wellbore trajectory model so asto generate trajectory data indicative of a trajectory of the wellbore;and determining, on the basis of the trajectory data, said one or moreoperating modes of the wellbore drilling system.
 11. The method of claim10, further comprising the steps of: automatically configuring acontroller of the wellbore drilling system with the one or moreoperating modes determined by the wellbore positioning system; andapplying the one or more operating modes.
 12. A wellbore positioningsystem arranged to determine the relative positions of a wellbore and anobject, the wellbore being represented by a first ellipse and the objectbeing represented by a second ellipse, wherein the first ellipserepresents a positional uncertainty of the wellbore and the secondellipse represents a positional uncertainty of the object, the systemcomprising: data receiving means arranged to receive input data relatingto a measured or estimated position of the wellbore and the object, theposition of the wellbore having a first set of parameters defining thefirst ellipse, and the position of the object having a second set ofparameters defining the second ellipse; expansion factor calculationmeans arranged to calculate an expansion factor representing an amountby which one, or both, of the first ellipse and the second ellipse canbe expanded with respect to one or both of respective first and secondsets of elliptical parameters so that the first and second ellipsesosculate, wherein calculating the expansion factor comprises determiningand solving a quartic equation that is based on the geometry of thefirst and second ellipses; and position determining means arranged todetermine, based on the calculated expansion factor, position dataindicative of the relative positions of the wellbore and the object. 13.The wellbore positioning system of claim 12, further comprisingoperating mode determining means arranged to determine, on the basis ofthe position data, one or more operating modes of a wellbore drillingsystem.
 14. The wellbore positioning system of claim 12, the systembeing operatively connected to a controller of the wellbore drillingsystem such that the controller of the wellbore drilling system isautomatically configured with the one or more operating modes determinedby the wellbore positioning system, the controller being arranged toapply the one or more operating modes.
 15. A computer program productcomprising a set of instructions which, when executed by a computingdevice, is configured to cause the computing device to carry out themethod according to claim
 1. 16. The computer program product of claim15, comprising a computer readable storage medium.